function [GG,RR,Cobs,eu,SDX,ZZ,initss,ssvec,flag,ssnames,namesstates,namesshocks,ACont,ALag,ALead,PSI] =...
    modExample3EQ(param,~,addsol)
% [GG,RR,Cobs,eu,SDX,ZZ,initss,ssvec,flag,ssnames,namesstates,namesshocks,ACont,ALag,ALead,PSI] =...
%    mod3EQ(param,~,addsol)
%
% This is an example mod file for usage with Alejandro's estimation code.
% I tried to pick the easiest model so that we could more easily understand
% the different sections of the code.  Further, there is only one shock and
% only one observable used, so some of the output may not make as much
% sense.  In particular, the deompositions by shock will be a little
% strange, since everything will be contributed to the one shock.
%
% Created 6/25/2015 by Todd Messer

% This is the basic 3 Equation NK DSGE model (Gali).  
%
% 1) Dynamic IS Equation:
%   \pi_t = \beta E_t\pi_{t+1} + \kappa y_t   
%
% 2) New Keynesian Phillips Curve
%   y_t = E_t y_{t+1} - \frac{1}{\sigma}(i_t - E_t\pi_{t+1})
%
% 3) Monetary Policy
%   i_t = \rho + \psi_{\pi}\pi_t + \psi_{y}y_t + \varepsilon^r_t
%
% 4) Natural Rate
%  	Set hte Natural RAte equal to 0 so there are no real shocks.
%
% 5) Monetary Policy Shocks
%   \varepsilon^r_t = \rho_r \varepsilon^r_{t-1} + \nu^r_t
%
% Preliminary Inputs and Outputs
% Here you can put in settings that will never change.  For example, this
% example model doesn't have any steady-state parameters to use in
% estimation, so we can set those to be empty.
initss  = []; 
ssvec   = []; 

% These flags may be old or not in use yet, so if you have any more
% specific questions about them ask Alejandro
flag.ssok   = 1;
flag.solver = 1;

%% Part 0: Modifications and Notes
% 6/25/2015 (Todd Messer): Added more comments and reorganized model code

%% Part 1: Parameters
% It is very important that you note the order of the model parameters.
% They will need to be inputted in the same order in your excel file that
% lists the  priors.  Typically, I would do the model parameters first,
% followed by the shock process parameters, followed by the standard
% deviations of the shocks.

% Model Parameters
betin   = param(1);
gamma   = param(2);
kappa   = param(3);
sigma   = param(4);
rho     = param(5);
psi_pi  = param(6);
psi_y   = param(7);

% Shock Parameters
rho_r   = param(8);

% Standard Deviations
STD_r   = param(9);

%% Part 2: Variables
% This is the ordering of the variables.  It is not important to note the
% order, but rather to make sure that each variable is assigned its own
% unique index.  I typically place the endogenous states first, followed by
% the shock processes, followed by the observables we need mapped.

% Endogenous States
y       = 1;
pi      = 2;
i       = 3;

% Shocks
shkr    = 4;

% Observables
yobs    = 5;

% Count - This can help to make sure you have every variable mapped into an
% index.  i.e. the largest index should match the numbers you put here.
NVarY   = 5;
Neq     = 5;

%% Part 3: Shocks
% The shock processes we put as endogenous variables correspond to, for 
% example, the AR processes that we have our shocks follow.  These rae the
% i.i.d. shocks that hit each period, and so we reset the index as they
% will be placed in the PSI matrix.
er      = 1;

% Count - Another check
NShocks = 1;

% SDX matrix - The standard deviation matrix.
SDX         = eye(NShocks);
SDX(er,er)  = STD_r;

%% Part 4: Model Declaration
% Here, we initialize the Lead, Lag, and Contemporeanous matrices for
% coefficients as well as the shock coefficient matrix.
ACont   = zeros(Neq,NVarY);    
ALead   = zeros(Neq,NVarY); 
ALag    = zeros(Neq,NVarY);  
PSI     = zeros(Neq,NShocks);  

%% Part 5: Steady-State 
% There may be some steady-state values that we want to calculate first.
% Reasons for doing so may be that they are used in multiple equations, or
% we want to output them at the end of the file to check.  An example may
% be the coefficient on marginal cost in the phillips curve - it may be
% helpful to see its magnitude.
bet         = (betin/100+1)^(-1);

%% Part 6: Equations
% Each index for a variable also correesponds to an equation.  That is, we
% need to have the same number of unknowns as equations.  Hence, when we
% assign the coefficients to each of the variables by equation, we can
% choose the row based on a variable.  At the end, each variable should
% have its own "equation."

% 1) IS Curve
% Describing this one, we can see that the 1st equation, or the "y"
% equation, has values in the contemporaneous matrix for output and the
% inflation, and values in the lead matrix for output.
ACont(y,pi) = 1;
ACont(y,y)  = -kappa;

ALead(y,pi) = -bet;

% 2) Phillips Curve
ACont(pi,y)  = 1;
ACont(pi,i)  = 1/sigma;

ALead(pi,pi) = -1/sigma;
ALead(pi,y)  = -1;


% 3) Monetary Policy
ACont(i,i)      = 1;
ACont(i,pi)     = -psi_pi;
ACont(i,y)      = -psi_y;
ACont(i,shkr)   = -1;

ALag(i,i)       = rho;

% 4) Shock
% Here, we can create the shock process.  This is an AR(1) process.
ACont(shkr,shkr)= 1;
ALag(shkr,shkr) = rho_r;
PSI(shkr,er)    = 1;

% 5) Measurement Y
% Sometimes the measurement equations don't map directly into one variable,
% but rather we have to create a new one.  For example, we feed in output
% growth as an observable, but output in the model is in log-linearized
% deviations fromt he mean.  Hence, we need to take the difference of the
% model output to put it in growth rates.
ACont(yobs,yobs) =  1;
ACont(yobs,y)    = -1;
ALag(yobs,y)     = -1;


%% Part 7: Solve Model
% We use the Anderson-Moore Algorithm.  Note that the vector of constants
% is just zeros as we have assumed a log-linearized model aroundmodmmmfdsam
% steady-state.  
[GG,~,RR,~,~,~,~,eu]=amasolve(ALead,ACont,ALag,zeros(NVarY,1),PSI); 

% Check for existence AND uniqueness.
flag.euok = isequal(eu(:),[1;1]);
if ~flag.euok
    ZZ=[];
    Cobs=[];
    SDX=[];
    RR=[]; 
    disp('Model has no solution')
    return
end

%% Part 8: Observables
% Here, we do an abbreviated version of the first part of the code.  That
% is, we assign indices to each of the observables, then we create a
% coefficient and output matrix to map the observables into some model
% variable.  
%
% In this model, we have only one observable and we map y_data into yobs. 

% Data Index
y_data  = 1;

% Count - As a check
NObs    = 1;

% Coefficient Matrix
ZZ = zeros(NObs,size(GG,2));
ZZ(y_data,yobs) = 1; 

% Constant Vector
Cobs        =   zeros(NVarY,1);
Cobs(yobs)  =   gamma;

%% Part 9: Names
% We input the names of each variable and shock.  It is important to note
% the order as the names of the states and shocks must be in the same order
% as how we inputted them originally.

% States
names.states  ={
    'Output'            % 1
    'Inflation'         % 2
    'Interst Rate'      % 3
    'MP Shock'          % 4
    'Output Growth'     % 5
    };

% Shocks
names.shocks={
    'Monetary Policy'   % 1
};

%% Part 10: Output
% Checks
% First, we check to make sure that we ahve the same number of state names
% as states themselves.  Then we do the same with the shocks.
if size(names.states)~=size(GG,1); 
    error('Wrong number of state names'); 
end 

if size(names.shocks)~=size(RR,2); 
    error('Wrong number of shock names'); 
end 

% Here's some actual output.  This model didn't have any steady-state
% values worth outputting so we leave those blank.
ssvec   = []; 
ssnames = {}; 

% Here's were we output the names.  I believe this could simplified
% earlier, but it may be easier to keep it this way.
namesstates         = names.states';
namesshocks.short   = names.shocks';
namesshocks.long    = names.shocks';
end